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ABSTRACT 

Strong, anisotropic turbulence reflecting magnetized filaments is considered, to model the diffusive 
acceleration of particles by shock waves in active galactic nucleus jets. We address that at knot A 
of the nearby M87 jet, the shock involving the filamentary turbulence can accelerate an iron nucleus 
to zetta-eV (ZeV; 10 21 eV) ranges. A smaller value of the particle diffusion coefficient is found to be 
essential to achieve a ZeV cosmic-ray accelerator, a "Zevatron." 

Subject headings: acceleration of particles — galaxies: individual (M87) — galaxies: jets — magnetic 
fields — methods: analytical — shock waves 



1. INTRODUCTION 

The morphology of self-organized filaments appears in 
a wide variety of celestial objects, particularly, in splen- 
did extragalactic jets (e.g., Cyg A: Perley et al. 1984; 
M87: Owen et al. 1989; 3C 273: Lobanov & Zensus 
2001), which extend to kiloparsec to megaparsec scales 
with narrow opening angles. The underlying physics 
seems to be ubiquitous; thereby, we encounter the sim- 
ilar, smaller scale structure in the Galactic center (GC; 
Yusef-Zadeh et al. 1984; Yusef-Zadeh & Morris 1987) 
and even in laboratory plasmas (Tatarakis et al. 2003). 
In particular, it is, at present, known that the GC 
region arranges more than 80 linear filaments (Yusef- 
Zadeh et al. 2004), including the nonthermal filaments, 
which probably embody the energetic flowing plasmas. 
Ordinarily, such filaments are permeated by, more or 
less, ordered magnetic fields, which participate in self- 
organizing the highly anisotropic structure. 

Relating to the issues of cosmic-ray acceleration, these 
observational results stimulate us to elaborate on a fun- 
damental theory of particle transport in the magnetized 
filaments, typical for astrophysical jets. As for the dif- 
fusive shock acceleration (DSA) of particles (e.g., Drury 
1983), the acceleration rate depends largely on the dif- 
fusion property, inferred from magnetic field strength, 
configuration, and turbulent state around the shock. Ac- 
cording to the conventional resonant scattering theory, 
the charged particles being bounded by the mean mag- 
netic field (or magnetic flux tube) are resonantly scat- 
tered by the fluctuating magnetic fields superposed on 
the mean field; and in this context, the DSA scenar- 
ios of extremely high energy (EHE) cosmic rays have 
been discussed for active galactic nucleus (AGN) jets, 
considering the parallel (Biermann & Strittmatter 1987; 
Rachen & Biermann 1993) and oblique (Honda & Honda 
2004) mean fields with respect to the shock normal di- 
rection. As the particles scatterer, the magnetohydro- 
dynamic (MHD) turbulence of Alfven waves with a Kol- 
mogorov spectrum was postulated, although the actual 
turbulent state in the jets has been unresolved to date. 
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At this juncture, recent polarization measurements us- 
ing very long baseline interferometry began to reveal that 
some AGN jets are permeated by the magnetic fields 
transverse to the jet axis (typically, 1803+784; Gabuzda 
1999). Over parsec scales, the fields are apparently so 
smooth that one cannot understand their polarity solely 
by invoking a trail of the fast-mode oblique shocks that 
refract the field lines to the off-axial direction. The fea- 
sible scenario responsible for the observed results is to 
allow huge currents launched from the central engine of 
the accretion disk (Honda & Honda 2002). As is theoret- 
ically known, the uniform currents are catastrophically 
unstable for the electromagnetic filamentation instability 
(Honda 2004 and references therein) that breaks up a 
single beam into many filaments carrying a limited cur- 
rent each (Honda 2000), and the nonlinear evolution re- 
sults in self-organizing toroidal magnetic fields (Honda 
et al. 2000a). In fact, large-scale toroidal magnetic fields 
have recently been discovered in the GC region (Novak 
et al. 2003), accompanied by many filaments. We ac- 
cordingly conjecture that the similar topology appears 
in extragalactic objects (Medvedev & Loeb 1999), inter 
alia, AGN jets. 

In this Letter, we put forth the filamentary AGN jets 
as a promising candidate for the cosmic-ray Zevatron 
(Blandford 2000). From a generalized quasi-linear trans- 
port equation coupled with the spectral intensity of the 
magnetized filamentary turbulence, we derive a diffusion 
coefficient of relativistic particles, and install it in the 
DSA model, in order to estimate the timescale of parti- 
cle acceleration. The remarkable feature is that, partic- 
ularly for high-Z nuclei, the coefficient becomes smaller 
than that derived from the resonant scattering theory in- 
voking the MHD turbulence, leading favorably to shorter 
acceleration time, viz., higher acceleration efficiency. We 
here demonstrate that in a feasible parameter domain for 
the brightest (radio to X-ray) knot, knot A (e.g., Wilson 
& Yang 2002), of the M87 jet, an iron nucleus can be 
accelerated to ZeV ranges. 

2. KINETIC THEORY OF PARTICLE DIFFUSION IN 
MAGNETIZED CURRENT FILAMENTS 

2.1. The Filamentary Jet Model 

Let us begin with Figure 1, illustrating a schematic 
view of our model for AGN jets as an EHE cosmic-ray 
accelerator. As now seems more likely, a jet, envisaged 
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Fig. 1. — Schematic of AGN jets, which operate as ZeV particle 
accelerators. The hierarchical structure of filaments illustrates that 
a jet, carrying huge current /, consists of numerous magnetized fil- 
aments with the smaller radial sizes (Honda & Honda 2002). In 
the transverse section of sample filaments, the "fish eyes" symbol- 
ically represent the hot currents pointing out of the page, and the 
loops with arrows the magnetic field lines. In the envelope region 
of bundled filaments, the large-scale toroidal field Bg can appear, 
establishing a magnetotail. 



as one or several filaments within the current resolution 
(Owen et al. 1989; Asada et al. 2000; Lobanov & Zcn- 
sus 2001), is hypothesized to be a bundle of numerous 
filaments, whose radial size is self-adjusted by an effec- 
tive Debye sheath (Honda et al. 2000a; Honda & Honda 
2002). For finite ion abundance, the charge conserva- 
tion law requires that the filamentary currents must carry 
the net charges corresponding to the number of electrons 
subtracted by that of positrons. In addition, there is di- 
rect observational evidence that the jet transports energy 
from the central engine to the radio lobe (Biretta et al. 
1995). Thus, the hot currents carried by the negatively 
charged electron-positron fluids flow in the opposite di- 
rection, as shown in Figure 1, while the cold return cur- 
rents compensating for positive ionic charges flow in the 
same direction via the filaments and/or external media 
(Honda & Honda 2002). 

From a macroscopic point of view, the vectors of the 
self-generated magnetic field appear to be randomly ori- 
ented in the transverse direction of the bundled filaments, 
such that their average is negligible, except for the enve- 
lope region. As far as the return currents flow, in part, 
via the external media (e.g., radio cocoons), the uncan- 
celled magnetic field pervades the exterior of the jet, 
establishing a long-range "magnetotail" of the toroidal 
field, in accordance with a pioneering observation that 
revealed not-so-weak (rather strong), ordered fields in 
the central kiloparsec-scale "hole" (cavity) of the inner 
radio lobe containing a jet (M87; Owen et al. 1990). 
The toroidal field will play a significant role in radially 
confining the jet (Honda & Honda 2002), as is consistent 
with the observational fact that pressure inside the jets 
is larger than that of the external media (Perley et al. 
1984; Owen et al. 1989). 



2.2. The Off-resonant Scattering of Particles by the 
Magnetized Filaments 

In a forest of the magnetized filaments, diffusion prop- 
erty of injected energetic particles is markedly different 
from that stemming from the small-angle resonant scat- 
tering (Drury 1983). When considering the DSA of the 
particles, it is necessary to know the effective diffusion 
coefficient for the direction normal to the shock front, 
K n . As outlined below, it can be derived from the Vlasov 
equation that describes the collisionless scattering of rel- 
ativistic test particles by the fluctuating magnetic fields 
B. For simplicity, the fields are supposed to be purely 
transverse to the z-axis, B = B_l(i", t), where r = (x,y), 
and electric fields are ignored, since they are preferably 
shorten out in the propagation time of jets (Honda et 
al. 2000a). Hereafter, the perpendicular and parallel 
subscripts refer to the z-axis reflecting the direction of 
the current filaments. We decompose the momentum 
distribution function of the particles into the averaged 
and fluctuating part, f p = (f p ) + Sf p , and combine the 
averaged Vlasov equation, (D/ p /D£) = 0, with the evo- 
lution equation of Sf p that contains the term involving 
(/p), linear to B^. For the manipulations, we carry out 
the Fourier transformation of all fluctuating quantities in 
terms of both r and t: T{r,t) = J d 2 kdwJ=" k , ClJ e i ( k r ^ wt ), 
where klz. As a result, we obtain 

1<£> = [ d 2 kd 2 k'd w d w v[( k+k ')--K-') t ] 

At c z \J 

x [v x (k' x A kV /)] 

•sl ^y •£}<*>)■ w 

where d/dt = d/dt + v • (d/dr), A kW ^w = A k (/) iW (oz 
denote the Fourier components of the vector potential 
that conforms to Bj_(r, t) = V x Am (r, t)z, c is the speed 
of light, and q and v are the charge and velocity of the 
particle, respectively. 

Even for the strong nonlinear fluctuations, the cor- 
relation function may be written as (A k)W A k ' )W ') — 
\A\^^S(k + k')S(uj + uj'), where the Dirac 5-function has 
been introduced (Tsytovich & ter Haar 1995). Assum- 
ing the stationary and homogeneous fluctuations, i.e., 
lo = —lu' and k = — k', and invoking the causality prin- 
ciple to handle the resonant dominator, we can express 
equation (1) in the form of the generalized quasi- linear 
equation. Especially for the typical current filamentation 
instabilities that include the Weibel instability, the un- 
stable mode is often quasi-static (Kazimura et al. 1998; 
Medvedev & Loeb 1999; Honda 2004) and isotropic on 
the transverse two-dimensional plane (Montgomery & 
Liu 1979), to give \A\i^ ~ 2\A\\5{uj), where k = |k|. 
Along these, we impose the off-resonance condition for 
wave-particle interaction, |k • v| 3> to, and take an aver- 
age over the pitch angle. Then, for example, the integral 
[cq. (1)] including the double partial derivatives of d/dp\\ 
reduces to 




where v±_ = \v±\ and = \p\\\. Here wc have defined 
the turbulent spectral intensity as Ik = 2nk (k 2 \A\\/ , 
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such that the magnetic energy density can be evaluated 
by u m = B 2 /(8tt) = J^Jkdk, where B 2 = (B 2 ± (t)). 
The minimum wavenumber fc m ; n is set to n/R, where R 
stands for the radius of a bundle of the filaments, namely, 
the "radius of the jet." 

The energy density of the fluctuations, w m , likely be- 
comes comparable to thermal pressure (Honda et al. 
2000a). In this sense, the anisotropic filaments may be 
regarded as the "strong turbulence." We see that the 
spectrum exhibits a power law, Ik oc fc~ Q , with its in- 
dex around a ~ 2 (Montgomery & Liu 1979), somewhat 
larger than «mhd ~ I — (5/3) for the classical MHD tur- 
bulence. The kinetic simulations also indicated that the 
filament coalescence led to the accumulation of larger fil- 
aments, releasing in part the free energy of flows (Honda 
et al. 2000b). These are compatible with the observed 
trends of the steepening of filamentary turbulent spectra 
(Carilli & Barthel 1996) and the merging of filaments 
(Owen et al. 1989). 

2.3. The Diffusion Coefficient 

For a simple estimation, let us assign the EHE particles 
having their energy of E — cp and the isotropic momen- 
tum distribution (f p ) oc oc E~@ with fi ~ 3 (Stecker 
& Salamon 1999; de Marco et al. 2003), so as to arrange 
equation (2) in the form of (d(/ p )/d£)|i = v e s,\\ (fp)i 
where v e ^ » reflects an effective collision frequency. The 
spatial diffusion coefficient can be then estimated as 
K|| ~ w|/(2i/ effi ||) ~ (V6ir/4)cE 2 /[/3(p-l)q 2 B 2 R], where 
we have used v 2 ~ v\/2 ~ c 2 /3 and the definition of u m . 
Furthermore, supposing the shock to propagate along the 
filaments (z-axis) gives n n = Kn. 

Here it may be instructive to look at the ratio of k„ 
to a conventionally used coefficient for omhd = 5/3, 
i.e., the Kolmogorov MHD turbulence (Bicrmann & 
Strittmatter 1987). For convenience, we compare the 
mean field strength of the turbulent magnetic fields 

in the MHD context to B [= ((B^)) 1/2 ] and intro- 
duce q = Ze, where Z and e are the charge num- 
ber of the particle and the elementary charge, respec- 
tively. Then we obtain the ratio of K n /«n,MHD ^ 
0.1[(£/1 ZeV)(26/Z)(l mG/B)(100 pc/R)} 5 ^. It is 
found that for kiloparsec-scale jets, larger values of Z 
and B lead to a smaller K n , thereby, favorably to a 
shorter DSA timcscale of ZeV particles (see § 3). Im- 
portantly, the gyroradius of injected energetic particles 
cannot be well defined in the present model, because of 
| (B_l) I < B in the interior of the jet (§ 2.1). When 
considering the particle confinement, therefore, we need 
to compare the three-dimensional rms deflection of the 
accelerated particle to the system size, rather than its 
gyroradius. Relating to this, we find kj_/k|| <~ 0(1), in 
contrast to a larger anisotropy that appears in the sim- 
ple MHD. It follows that the radial size, which is smaller 
than the length of jet, affects the confinement. Note here 
that the radially decaying magnetotail may play an addi- 
tional role in confining the leaky energetic particles with 
their long mean free path (mfp) of Aj_ (~ v± / v c s,±) ~ R- 

3. ZEV ACCELERATION OF THE IRON NUCLEUS IN THE 

M87 JET 

We are concerned with nearby radio galaxy M87 
(Virgo A) residing in the center of the Virgo Cluster, 



which is known as a sub-TeV (possibly TeV) gamma-ray 
emitter (Aharonian et al. 2003). For application to the 
jet in the core of the giant elliptical galaxy, one issue 
worth noting is that a condensation of heavy elements 
was discovered in the core (Gastaldcllo & Molendi 2002). 
Presuming that the abundance of heavy elements includ- 
ing iron is finite in the jet, as confirmed in a Galactic 
microquasar jet (Kotani et al. 1996), as well as referring 
to up-to-date observational results around ~ 100 EeV 
(10 20 eV), compatible with an assumption of heavier 
composition primaries (Ave et al. 2000; Risse et al. 
2004), we pay attention to the acceleration of iron nu- 
clei (and protons). 

Taking account of the Fermi type I mechanism that is 
plausible for the knotlike regions of the M87 jet (Heinz 
& Begelman 1997), the mean acceleration time can be 
expressed as i acc ~ [3/(Ui - U 2 )][{k dA /U 1 ) + (k d ^/U 2 )}, 
where the subscripts i = 1,2 indicate the upstream and 
downstream region of the shock and Ui are the flow 
speeds in the shock rest frame (Drury 1983). When as- 
suming the values of a, (3, and R are constants, we have 
K n ,i = K n .2, because of the relation of B\ = B 2 derived 
from the condition that the current density Jnz must be 
continuous across the shock front. Then the above ex- 
pression of the acceleration time reduces, in the strong 
shock limit (Biretta et al. 1991), to i acc ~ 20« n /[/ 2 , 
where K n = K n .i, and U = U\ is the shock speed in the 
laboratory frame. 

Concerning the energy constraints relevant to the tem- 
poral scale, the shortest timescale that most severely re- 
stricts the acceleration of nuclei is arguably the shock 
propagation time, which is estimated as i s h ~ L/U, 
where L represents the propagation distance of the shock. 
It is noted that for the kiloparsec-scale jet, the syn- 
chrotron loss of the accelerated nuclei is ignorable for 
the inherent weak S-field, and radial adiabatic expansion 
as well is ineffective for the self-collimating jet (Honda 
& Honda 2002, 2004). For the moment, balancing i acc 
with i s h is adequate for estimating the achievable max- 
imum energy, and then, solving for E yields the scaling 
of E m , t oc ZB(LRU) 1 / 2 (for the temporal limit). In ad- 
dition, the particle confinement radius, i? c , limits the 
acceleration in terms of the spatial scale. As remarked 
in § 2.3, we expect R c > R, taking the effects of the 
magnetotail into consideration. By simply equating the 
transverse mfp Aj^ with R c , we obtain another scaling of 
the achievable maximum energy, £ miS oc ZBp 1 l 2 R (for 
the spatial limit), where p = R c /R. As a consequence, 
in the optically thin regions distant from the central en- 
gine, the actual maximum energy of accelerated nuclei 
can be expressed as E m = min(.E miS , E m ^). Note that 
the larger R leads to the higher E m , because the accu- 
mulating larger filaments more largely contribute to the 
deflection of particles. This R dependence is quite con- 
trary to that deduced from the MHD model, in which 
a larger resonant gyroradius involves a longer coherence 
time, that is, longer acceleration time (Honda & Honda 
2004). 

For the evaluation of E m , t , we unfix the value of U, 
since it seems to be uncertain yet, although the proper 
motion of some knots was surveyed in detail and found to 
exhibit average speeds around 0.5c (Biretta et al. 1995). 
An upper limit of L is estimated by using an ad hoc rela- 
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Fig. 2. — Maximum energy of the accelerated iron nucleus 
E m \z=2& vs. shock speed U in units of c, for the M87 jet/knot A. 
The rms strength of the random magnetic fields. B, is set to (a, 
6) 4.59 mG (Owen et al. 1989), (c) 640 fj,G (Meisenheimer et al. 
1996), (d) 310 fiG (Owen et al. 1989), and (e) 28 ;uG (Heinz & 
Begelman 1997). The viewing angle 6 and the radial extension 
factor p are chosen as (a) 30° (Bicknell & Begelman 1996) and 
2, or (6— e) 42.5° (Biretta et al. 1995) and 1, respectively. Here 
the opening angle (f> and the radius R have been fixed to 6.9° 
(Reid et al. 1989) and 55 pc (Biretta et al. 1991), respectively. 
The hatched area indicates E m in a plausible parameter domain 
[between (6) and (d)]. For comparison, we also indicate the pattern 
speeds of the bright edge of the knot (0.42 — 0.50; a band with 
a fat bar, where the darkly shaded zone suggests the most plau- 
sible E m ) and the entire knot (> 0.454; arrow; Biretta et al. 1995). 

tion, L/R ~ 360 /ircj) sin 9, for the opening angle (j) ss 6.9° 
(Reid et al. 1989) and viewing angles 6 sa 42.5° (Biretta 
et al. 1995) or 30° (e.g., Bicknell & Begelman 1996) of 
the M87 jet. As for E m>a , the value of p can presumably 
take ~ 10 at most (Owen et al. 1990). In Figure 2, for 
R = 55 pc at the brightest knot, knot A (Biretta et al. 
1991), given p, L(0), and B, now we plot E m \z=26 as a 
function of U in units of the speed of light c. It is found 
that even for a conservative value of p — 1 (Figs. 2b— 2e), 
the particle acceleration is, in the region of U < 0.42 for 
L(42.5°) « 1.4 kpc, not limited spatially, but temporally, 
whereas in the region of U > 0.42, vice versa. For a con- 
vincing value of B = 640 pG (Fig. 2 c; Meisenheimer et 
al. 1996), which is in a plausible range of B — 310 pG to 
4.59 mG (Figs. 2b-2d; hatched area; Owen et al. 1989), 
E m \z=26 ~ 1 ZeV is fairly achieved for U > 0.42, al- 
though in a marginal value of B = 28 pG (Fig. 2e; Heinz 
& Begelman 1997), it reduces to ~ 50 EeV. In an op- 
timistic parameter set of p = 2, L(30°) ~ 1.8 kpc, and 



B = 4.59 mG (Fig. 2a), E m \ z =26 ~ 10 ZeV is achievable 
for U ~ 0.5, whereby E m \ z =i ~ 400 EeV. 

The validity of this DSA model could be corroborated 
by the complementary calculation of the maximum en- 
ergy of an accelerated electron. For an electron with en- 
ergy much lower than E m , an eminent increase of acceler- 
ation efficiency is expected, on account of K n /K n ,MHD <C 
1. In regard to the energy limit, it is sufficient to consider 
the synchrotron loss (for M87; Biermann & Strittmatter 
1987) with its timescale of i syn oc 7~ 1 _B~ 2 , where 7 is 
the Lorentz factor of the electron. Taking the balance of 
iacc for Z = \ with i syn yields the expression of the max- 
imum Lorentz factor, which scales as 7 m oc R^U 2 / 3 . 
A noticeable thing here is that both i acc and t syn have 
a common dependence of being proportional to B~ 2 , so 
that B dependence of 7 m disappears, in contrast to the 
expression proposed by Biermann & Strittmatter (1987). 
The resultant maximum Lorentz factor 7 m ~ 10 11 , pre- 
dicted for R = 55 pc and U ~ 0.5, is amenable to the 
recent observational results that indicate no steep syn- 
chrotron cutoff even in the X-ray band of 10 17 — 10 18 Hz 
(reflecting 7 - 10 7 - 10 s : Marshall et al. 2002; Wilson 
& Yang 2002), considerably higher than the previously 
suggested cutoffs of ~ 10 15 Hz (e.g., Meisenheimer et al. 
1996; Heinz & Begelman 1997). 



4. CONCLUSIONS 

We have constructed a scenario of strong, off-resonant 
scattering of test particles by the quasi-static magne- 
tized filaments created via the current filamentation 
instability, to model the DSA in AGN jets. For the M87 
jet as an example, considering the most severe energy 
restriction from the spatiotemporal scale, we have 
estimated the achievable highest energy of accelerated 
particles at knot A. The results indicate that there is 
a wide range of plausible parameters where a shock 
can energize an iron nucleus to ~ 1 ZeV (whereby a 
proton to ~ 40 EeV) and more. At the moment, the 
problem of intergalactic transport of particles remains 
unsolved, though M87 is nearby (about 16 Mpc from 
us), so that the effects of a collision with the microwave 
background are not so significant (Stecker & Salamon 
1999). The present consequences might provide a key 
to elucidate the origin of EHE cosmic rays with energy 
above 100 EeV, in the context of the point-source model. 
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